float my_rand(unsigned int *seed) {
	unsigned long a = 16807; // constants for random number generator
	unsigned long m = 2147483647;
	// 2^31 - 1
	unsigned long x = (unsigned long) *seed;
	x = (a * x)%m;
	*seed = (unsigned int) x;
	return ((float)x)/m;
}

__kernel
void calcPI (__global float *gpu_pi) {

	int TRIALS_PER_THREAD=4096;

	unsigned int tid = get_global_id(0);
	int points_in_circle = 0;
	float x, y;
	unsigned int seed = tid + 1;

	// starting number in random sequence
	int i;
	for(i = 0; i < TRIALS_PER_THREAD; i++) {
		x = my_rand(&seed);
		y = my_rand(&seed);
		points_in_circle += (x*x + y*y <= 1.0); // count if x & y is in the circle.
	}
	gpu_pi[tid] = 4.0 * points_in_circle / (float) TRIALS_PER_THREAD;   
}
